On the qualitative study of a discrete-time phytoplankton-zooplankton model under the effects of external toxicity in phytoplankton population

The current manuscript studies a discrete-time phytoplankton-zooplankton model with Holling type-II response. The original model is modified by considering the condition that the phytoplankton population is getting infected with an external toxic substance. To obtain the discrete counterpart from a continuous-time system, Euler's forward method is applied. Moreover, a consistent discrete-time phytoplankton-zooplankton model is obtained by using a nonstandard difference scheme. The boundedness character for every positive solution is discussed, and the local stability of obtained system about each of its fixed points is discussed. The existence of period-doubling bifurcation at a positive equilibrium point is discussed for the discrete system obtained by Euler's forward method. In addition, the comparison of the consistent discrete-time version with its inconsistent counterpart is provided. It is proved that the discrete-time system obtained by using a nonstandard scheme is dynamically consistent as there is no chance for the existence of period-doubling bifurcation in that system. In order to control the period-doubling bifurcation and Neimark-Sacker bifurcation, an improved hybrid control strategy is applied. Finally, we have provided some interesting numerical examples to explain our theoretical results.


Introduction
Many scientific investigations have focused on mathematical modelling and studying population dynamics over the last two decades. Moreover, in that period, some well-known mathematical models in abstract ecology have been presented, and the most famous one is the Lotka-Volterra model [1]. The study of creature movement and scattering has turned out to be a vital element for understanding a sequence of environmental interrogations connected with the spatiotemporal analysis of population dynamics [2]. Plankton is enormously flexible in abundance, both temporally and spatially. Plankton discrepancy depends on the regular and the spatial structure's physical procedure. Biological processes include, for instance, growth, eating and behaviour. Moreover, physical processes are involved as adjacent stirring. The nonlinearity of bio networks contributes entirely to the spatial union in plankton distributions [3]. In marine ecology, the word plankton refers to a class of freely moving and hazily swimming organisms. Generally, plankton is of two kinds: phytoplankton and zooplankton. Phytoplankton classes have a petite size with a single-celled structure [4]. By primary formation, sinking and decease, they put forth a worldwide effect on the environment by efficiently carrying 2 from the ocean's exterior layer to the deep aquatic sediments. The algal species grow in large quantities in marine and limnic atmospheres. The stages of enhanced development slowed stagnation and rapid degeneration of cell counts collectively creating an algal bloom. This type of rapid variation in the phytoplankton population density is a uniqueness in the plankton ecosystem [4]. Even though the unexpected appearing and fading of blooms are not well understood, the opposing properties of damaging algal blooms on the health of humanity, marine population, fisheries industry, and tourism are proven (see [4]). The analysis of the interaction of phytoplankton and zooplankton on the happening of bloom is of interest to several scientific studies [5]. Phytoplankton creates venoms to evade the predation by zooplankton. Decrease of eating force of zooplankton in arrears to the ejection of venomous materials by phytoplankton is one of the motivating topics of investigation in the past years. The venom making phytoplankton decreases the grazing force on them. Thoughtful of the dynamics of plankton inhabitants and their relations are of primary significance. Numerous phytoplankton species are poisonous to zooplankton. This poisonousness affects the spreading of phytoplankton and zooplankton populations. Scientific modelling of plankton population is an adequate substitute procedure in enlightening our information about the biological and physical processes connecting to the plankton ecosystem [5]. In [6], the authors have discussed a plankton-nutrient system related to marine ecology by examining planktonic blooms. The authors in [7] have discussed the effects of seasonality on the plankton population and explored the impact of periodicity on their dynamical behaviour.
By taking the intense depiction of viral diseased phytoplankton and infections, the author, in [8] has developed and studied two different mathematical models for the plankton population. The co-occurrence of competitive predators and possible extinction of prey due to the effects of predation are discussed in [9]. Above and beyond, they have studied a two-zooplankton one-phytoplankton model together through the observation of harvesting.
A nutrient-phytoplankton model is investigated in [10] to scrutinize the dynamics of phytoplankton blooms. Yunfei et al. [11] have presented a biological model for zooplankton-phytoplankton dynamics by including an extra condition of harvesting of populations. In addition, they have explained that over-manipulation may act as a source of the destruction of the population; however, proper harvesting must confirm the defensibility of the population. Similarly, several researchers have studied the dynamics of phytoplankton-zooplankton interacting systems by using the cooccurrence of plankton, a nutrient cause, the influence of harvesting, or the vicious influence of the plankton system [9,10,11,12,13,14,15]. Chattopadhyay and Sarkar [16] have presented and studied a mathematical model by considering the time delay in toxin escape by phytoplankton. In addition, it is suitable to lead with the poison-producing delay throughout the learning of the dynamical behaviour of phytoplankton-zooplankton population models. The efforts made in [17,18,19,20,21,22,23] inspired us to study and discuss the dynamical behaviour of a phytoplanktonzooplankton inhabitant's model using poisonousness. Additionally, this poisonous material is released by phytoplankton and sometimes by an external toxic source. For the study of some interesting models in mathematical ecology and mathematical chemistry, we refer a reader to [24,25,26,27]. Here, we have considered the basic idea of modelling phytoplankton-zooplankton interactions, which is taken from Zhang and Rehim [28]. Moreover, while studying our mathematical model, we focus on the following situations.
• We consider that ( ) and ( ) are densities of zooplankton and phytoplankton populations at any time .
• The zooplankton population is continuously consuming the phytoplankton population and recycles them into their community.
• We have considered that the zooplankton population becomes diseased by ingesting the diseased phytoplankton population. Furthermore, the contagion in phytoplankton may produce due to some external toxic substances (see [28]). • The phytoplankton has exponential growth in the nonappearance of the zooplankton population. Where is their logistic rate of change, and is the maximum carrying capacity of the surroundings [29]. • We have taken that the time lag for phytoplankton's creation and mediation of toxic stuff is zero.
Under these conditions, as mentioned earlier, we have the following phytoplankton-zooplankton model from [28]: Kaung [30] examined the behaviour of a prey-predator model by using the Holling type-II response function. Moreover, he explained that studying the dynamical assets of prey-predator models using Holling type response is preferable to studying dynamics of predator-prey models without using Holling response [31]. Generally, we define the Holling type-II response function by the basic definition of rectangular hyperbola, and its mathematical form is: Where is some constant, by using Holling type-II response [31], we have the next mathematical shape of the system (1): • The functional response ( ) ( ) + ( ) characterizes the rate at which phytoplankton population is eaten by zooplankton. In addition, it causes a rise in the development rate of zooplankton and this development rate is represented by ( ) ( ) + ( ) . • The term 2 ( ) in the system (2) shows that the infection developed in the phytoplankton population is due to external toxic substances.
Additionally, the term 2 2 ( 2 ) = 2 > 0 illustrates the fast-tracking growth of poisonous substance corresponding to phytoplankton population, as approximately every type of phytoplankton population is increasingly consuming toxic substances. Where the parametric values in the system (2) are non-negative and defined along these lines: • : the highest seizure rate of zooplankton on phytoplankton.
• : the poisonousness rate of phytoplankton for every unit biomass. • : the decease rate of zooplankton inhabitants. According to Strogatz [32], chaos occurs in a continuous system when it is at least 3-dimensional. Therefore, it is clear that chaos ceased to exist in the system (2). However, in the case of a discrete-time map, chaos can be observed in 1-dimension. Hence, there is a chance for the existence of chaos in the 2-dimensional discrete-time phytoplankton-zooplankton model. In addition, if chaos does not exist after the discretization then the discretized system is said to be dynamically consistent. Motivated by the aforementioned prosperous characteristics of discrete-time dynamical systems, it is necessary and exciting to study the qualitative behaviour of the discrete-time version of the system (2). Moreover, our main aim is to study a consistent counterpart of the system (2) such that there is a minimal change in the dynamical behaviour of the discretized system compared to the original continuous system. Therefore, by using Euler's forward method with step size , we have the following discrete-time version of (2): Moreover, to obtain a consistent counterpart of system (3) and by applying the Micken's type nonstandard scheme on the model (3), we get the following discrete-time mathematical model (see [33]): The next part of this manuscript is structured along these lines: • The boundedness character for every positive solution of system (4) is discussed in section 2.
• The existence of fixed points and local stability of system (3) and (4) about each of them is investigated in section 3.
• The presence of Neimark-Sacker bifurcation about the unique positive fixed point of system (4) is discussed in section 4.
• The existence of period-doubling bifurcation about one and only positive fixed point of system (3) is discussed in section 5.
• We discussed a modified hybrid control strategy for controlling the chaos, period-doubling bifurcation and Neimark-Sacker bifurcation in section 6. • A comprehensive numerical simulation is provided in section 7 to support each theoretical investigation.

Boundedness of solutions of system (4)
We assume that 0 > 0 and 0 > 0, with ≥ for all ≥ 0. Then every solution ( , ) of the system (4) must fulfils > 0 and > 0 for all ≥ 0. Moreover, from first part of (4) we get On solving (5), and by applying limit we get for every ≥ 0. In addition, for ≥ for all ≥ 0, from second part of system (4), we get ) .
Hence, one can obtain the upper bound for zooplankton population as for all ≥ 0. Finally, we have the following theorem about the boundedness of all solutions of (4).

Lemma 3.1. [34] Let
be the characteristic equation obtained from a 2 × 2 jacobian matrix . Moreover, be the jacobian matrix of system (3) or (4) about each of its equilibrium point. Additionally, assume that and are respectively trace and determinant of and (1) > 0. Then: As 1 and 2 are eigenvalues of (9) then we have the subsequent topological type outcomes related to the stability of . Where p be any arbitrary fixed point of (9). The point is identified as sink if | 1 | < 1 and | 2 | < 1, it is locally asymptotically stable. The point is known as source if | 1 | > 1 and | 2 | > 1, as . The point is acknowledged as non-hyperbolic if either 4 or 5 is satisfied.
Next, by considering the fixed point of the system (4) we get the following result: be the jacobian matrix of system (4) about ( + , 0). Then, the fixed point ( + , 0) of system (4) will be classified as In Fig. 1 the topological classification of fixed point ( + , 0) for some values of parameters is given. Now, we assume that (8) remains true, then by considering the fixed point of system (4) we have Moreover, the characteristic polynomial obtained from matrix [ ] is given as follows: where, ) .
By assuming that (8) remains true and by performing some mathematical operations, it follows that: and Remark 3.1. Assume that (8) remains true, then there is no chance of period-doubling bifurcation in system (4) as (−1) > 0 for every , , , , , , > 0.
We have the following theorem for the possible validation of Remark 3.1. Proof. Assume (8) and let be the positive fixed point of system (4) and Then, (−1) > 0 if and only if we have As then, we have Formally, from (15) we get The inequality (16) is true for every choice of , , , , , and ≥ . Which completes the proof. □ The following proposition around the local stability analysis of the system (4) about .

The point undergoes the Neimark-Sacker bifurcation if and only if
In Fig. 2, the topological classification of fixed point for some parametric values is given. Moreover, it can be seen that increasing the parameter may cause an increase in the bifurcation region.

Neimark-Sacker bifurcation for system (4)
In this section, we study the existence of the Neimark-Sacker bifurcation about the only positive steady-state of the system (4). We have used the standard theory of bifurcation for the direction and presence of this kind of bifurcation. In recent times, Neimark-Sacker bifurcation associated with some discrete-time mathematical systems has been examined by several scientists [35,36,37,38]. Additionally, when mathematical models are in differential form, we refer to [39,40,41] for specific considerations associated with Hopf bifurcation. First, we confirm that the positive equilibrium point of the system (4) experiences the Neimark-Sacker bifurcation whenever is chosen as a bifurcation parameter. One can see from Proposition 3.3 that roots 1 , 2 of (8) are complex and satisfy Furthermore, under the supposition that (8) remains true, we have the following set: Then, the positive fixed point of system (4) experiences the Neimark-Sacker bifurcation for parameter when it varies in a slight neighbourhood of ̂, which is given aŝ In addition, assume that ( , , , , , , , ) ∈ Φ 2 then the system (4) is characterized equivalently with the following planer map: To discuss and analyze the normal form theory for Neimark-Sacker bifurcation for fixed point , of (28), we suppose that 1 represents a small perturbation in ̂. Then the perturbed mapping for (19) can be described by the next map: where,  The characteristic equation for the jacobian matrix of system (21) calculated at (0, 0) is given as: where,̂( ) .
Hence, the mapping centred to set ℧ (0, 0, 0) is defined in such a way: where, Next, we have the following real numbers: ) Hence, by the study as mentioned earlier, we have the following result connected to the presence and direction of period-doubling bifurcation for mathematical system (3) about ̄.
Theorem 5.1. Assume that the parameter changes in least neighbourhood of ̄and 12 ≠ 0, then the system (3) experiences the period-doubling bifurcation at one and only positive equilibrium ̄. Additionally, the period-two trajectories that bifurcates from ̄are stable for 12 > 0, and if 12 < 0, then these trajectories are unstable.

Chaos control
To control inconsistent, accidental and irregular behaviour in any biological system, chaos control is well thought-out to be a practical tool for avoiding this complex and disordered behaviour. For additional details associated with the biological significance of chaos control and its applied use in the real world, we mention to [45]. In this part of the manuscript, we use a simple chaos control method for the system (3). Furthermore, there are many chaos control techniques for discrete dynamical systems. We refer a reader to [43,44,45,46,47,48,49] for additional details connected to these methods. We implement a generalized hybrid control technique (see [33,34,35,36,37]). The generalized hybrid control method [48] is centred on parameter perturbation and a state feedback control technique. By implementing generalized hybrid control methodology (with control parameter Θ ∈ (0, 1]) to the system (3), we get: Then, system (34) is controllable provided that its fixed point , is locally asymptotically stable. Additionally, the jacobian matrix for system (34) at its positive fixed point is calculated as follows: Formerly, the characteristic polynomial for the above-mentioned jacobian matrix is specified by  . Then, the mathematical system (4) experiences the Neimark-Sacker bifurcation whenever the bifurcation parameter certainly passes through = 0.4115075. Furthermore, for these values, bifurcation diagrams and MLE are shown in Fig. 3 (a)-(c). In addition, some phase portraits for system (4) are shown in Fig. 4  ] .
Additionally, the characteristic equation from 3 is given as The complex roots for (35)      ] .
Additionally, the characteristic equation from 3 is given as The roots for (36) are 1 = −1, and 2 = 199.87637628645194. Moreover, the bifurcation diagram for and MLE are shown in Fig. 9 (a)-(b).

Conclusions
In our work, we have briefly explained the dynamics of a continuous-time phytoplankton-zooplankton model [28]. Initially, the model taken from [28] is modified by considering the situation that some outer toxic substances are infecting the phytoplankton population. Furthermore, these poisonous substances have an accelerating progression compared to the phytoplankton inhabitants. Then, we have obtained a discrete-time type of model (2) by using Euler's forward method. According to Strogatz [32], chaos occurs in a continuous system when it is at least 3-dimensional. Therefore, it is clear that chaos ceased to exist in the system (2). However, in the case of a discrete-time map, chaos can be observed in 1-dimension. Hence, there is a chance for the existence of chaos in the 2-dimensional discrete-time phytoplankton-zooplankton model. In addition, if chaos does not exist after the discretization then the discretized system is said to be dynamically consistent. Moreover, our main aim is to study a consistent counterpart of the system (2) such that there is a minimal change in the dynamical behaviour of the discretized system compared to the original continuous system. Hence, using the nonstandard finite difference scheme, we have obtained a consistent counterpart (4) of the system (2). The boundedness of every positive solution (4) is discussed. The local asymptotic stability of obtained mathematical system (4) is discussed about each of its fixed points. Moreover, to show the complex behaviour in the mathematical system (4), the presence of Neimark-Sacker bifurcation about a positive fixed point is discussed (see Fig. 3 (a)-(c)). It is shown that the system (4) is dynamically consistent, and there is no chance for the existence of period-doubling bifurcation for system (4) (see Theorem 3.1). The existence of period-doubling bifurcation for a unique positive fixed point of the system (3) is shown mathematically, and some exceptional numerical examples are provided (see Fig. 9 (a)-(b)). It is demonstrated numerically that the system (3) experiences the Neimark-Sacker bifurcation (see Fig. 7 (a)-(c)) and period-doubling bifurcation for an extensive range of stepsize . From Fig. 8 (e) and Fig. 8 (f), the existence of chaos in the system (3) can be seen, which proves the inconsistency of Euler's method (see Fig. 8 (e)-(f)). Finally, Neimark-Sacker bifurcation and period-doubling bifurcation are effectively controlled by using a generalized hybrid method. It is shown that when the generalized hybrid method is applied to the system (3), it regains the system's stability (3) for the maximum range of control parameters (see Fig. 10 (a)-(b)). In addition, the generalized technique of control is comparatively more effective than the hybrid method [49]. It is centred on response control and brings back the system's stability for an extensive range of parameters (see Fig. 11). Moreover, from a numerical study, it is realized that the modified hybrid control methodology is appropriate for controlling the Neimark-Sacker bifurcation and chaos.

Author contribution statement
Muhammad Salman Khan: Analyzed and interpreted the data; Wrote the paper. Maria Samreen: Contributed reagents, materials, analysis tools or data. J.F. Gomez-Aguilar: Analyzed and interpreted the data; Wrote the paper. Eduardo Perez-Careta: Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement
This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Data availability statement
No data was used for the research described in the article.

Declaration of interests statement
The authors declare no competing interests.

Additional information
No additional information is available for this paper.